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Abstract 

In this paper we revisit the Brownian motion on the basis of the fractional Langevin equation 
which turns out to be a particular case of the generalized Langevin equation introduced 
by Kubo in 1966. The importance of our approach is to model the Brownian motion more 
realistically than the usual one based on the classical Langevin equation, in that it takes into 
account also the retarding effects due to hydrodynamic back-flow, i.e. the added mass and 
the Basset memory drag. We provide the analytical expressions of the correlation functions 
(both for the random force and the particle velocity) and of the mean squared particle 
displacement. The random force has been shown to be represented by a superposition of 
the usual white noise with a "fractional" noise. The velocity correlation function is no 
longer expressed by a simple exponential but exhibits a slower decay, proportional to t~ 3 / 2 
for long times, which indeed is more realistic. Finally, the mean squared displacement is 
shown to maintain, for sufficiently long times, the linear behaviour which is typical of normal 
diffusion, with the same diffusion coefficient of the classical case. However, the Basset history 
force induces a retarding effect in the establishing of the linear behaviour, which in some 
cases could appear as a manifestation of anomalous diffusion to be correctly interpreted in 
experimental measurements. 



1 Paper published in Modern Problems of Statistical Physics, Vol 8, pp. 3-23 (2009): a journal founded to 
the memory of Prof. Ascold N. Malakhov, see http://www.mptalaj11.0rg/i.html 
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1 Introduction 

Since the pioneering computer experiments by Alder k Wainwright (1970) [2], which have 
shown that the velocity autocorrelation function for a Brownian particle in a dense fluid 
goes asymptotically as t _3//2 instead of exponentially as predicted by stochastic theory, many 
attempts have been made to reproduce this result by purely theoretical arguments. 

Usually, hydrodynamic models are adopted to generalize Stokes' law for the frictional force 
and obtain a Generalized Langevin Equation (GLE)), see e.g. Zwanzig k Bixon (1970, 1975) 
[H511HS], Widom (1971) [S3], Case (1971) [12], Mazo (1971) [53], Ailwadi k Berne (1971) PQ, 
Nelkin (1972) [60], Hynes (1972) [35], Chow k Hermans (1972a, 1972b, 1972c) 031 EH Q3] , 
Hauge k Martin-L6f (1973) [35], Dufty (1974) [21], Bedeaux k Mazur (1974) [5], Hinch 
(1975) [37], Pomeau k Resibois (1975) [65], Warner (1979) [82], Reichl (1981) [68J, Paul k 
Pusey (1981) [62], Felderhof (1991) [27], Clercx k Schram (1992) [IS]. 

Recently, a great interest on the subject matter has been raised because of the possible 
connection among long-time correlation effects, fractional Langevin equation and anomalous 
diffusion, see e.g. Muralidar et al. (1990) [59], Wang k Lung (1990) [80], Wang (1992) [79], 
Porra et al. (1996) [66], Wang k Tokuyama (1999) [ST], Kobelev k Romanov [42], Metzler 
k Klafter [55], Lutz (2001) [45], Bazzani et al. (2003) [4], Budini k Caceres (2004) [7J, 
Pottier k Mauger [67], Coffey et al. (2004) [UJ, Khan k Reynolds (2005) [39], Fa (2006) 
[24], Vinales k Desposito [76], Fa (2007) [25], Vinales k Desposito [77], Taloni k Lomholt 
(2008) [73] , Vifiales et al. [78], Desposito k Vinales (2009) [IE], Figueiredo Camargo et al. 
(2009a, 2009b) [283 I2H] , Eab & Lim (2009 ) [22]. 

We recall that anomalous diffusion is the phenomenon, usually met in disordered or 
fractal media see e.g. Giona k Roman (1992) [31], according to which the mean squared 
displacement (the variance) is no longer linear in time (as in the standard diffusion) but 
proportional to a power a of time with < a < 1 (slow diffusion) or 1 < a < 2 
(fast diffusion), see e.g. Bouchaud k Georges (1990) [6]. In general anomalous diffusion 
phenomena are related to generalized diffusion equations containing fractional derivatives 
in space and/or in time and, on this respect, the literature is huge. We limit to recall 
some (review) papers and books that have mainly attracted our attention in the last decade 
to which the reader can refer along with the references therein: Carpinteri k Mainardi - 
Editors (1997) [UJ, Saichev k Zaslvasky (1997) [71], Grigolini et al. (1999) [34], Uchaikin 
k Zolotarev (1999) [75], Metzler k Klafter (2000) [56], Hilfer - Editor (2000) [36], Mainardi 
et al. (2001) 09] Gorenflo et al. (2002) [33], Saichev k Utkin (2002) [70], Meerschaert et al. 
[54], Metzler k Klafter (2004) [57], Piryatinska et al. (2005) [63], Dubkov et al. (2008) [22], 
Uchaikin (2008) [74], Klages et al. - Editors (2008) [4T]. 

We also point out that Kubo (1966) [43] introduced a generalized Langevin equation (GLE), 
where the friction force appears retarded or frequency dependent through an indefinite 
memory function. To be consistent with the fluctuation- dissipation theorem, in GLE the 
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random force is no longer a white noise (as in the classical Langevin equation) with a 
frequency spectrum related to that of the velocity fluctuations. As a matter of fact, the 
hydrodynamic models appear as particular cases of GLE, as noted by Kubo et al. (1991) 

m- 

In this paper we revisit the Brownian motion via a generalized Langevin equation of fractional 
order based on the approach started in earlier 1990 's by Mainardi and collaborators, see 
Mainardi et al. (1995) [51] and Mainardi & Pironi (1996) |50j. We model the Brownian 
motion by taking into account not only the Stokes viscous drag (as in the classical model) 
but also the retarding effects due to hydrodynamic back-flow, i.e. the added mass, and the 
Basset history force. 

The plan of the paper is the following. After having reviewed in Section 2 the classical 
Brownian motion, in Section 3 we extend the theory according to the hydrodynamic 
approach. On the basis of the fluctuation-dissipation theorem (see Appendix A) and of the 
techniques of the fractional calculus (see Appendix B), we provide the analytical expressions 
of the autocorrelation functions (both for the random force and the particle velocity) and 
of the displacement variance. Consequently, the classical theory of the Brownian motion 
is properly generalized. Significant results are shown and discussed in Section 4, where we 
point out different diffusion regimes. Finally, conclusions are drawn in Section 5. 



2 The classical approach to the Brownian motion 

We assume that the Brownian particle of mass m executes a random motion in one dimension 
with velocity V = V(t) and displacement X = X(t). The classical approach to the Brownian 
motion is based on the following stochastic differential equation (Langevin equation), see e.g. 
Kubo et al. (1991) [S], 

dV 

m — = F v {t) + R{t) , (2.1a) 

where F v (t) denotes the frictional force exerted from the fluid on the particle and R(t) 
denotes the random force arising from rapid thermal fluctuations. R(t) is defined such that 
for any t\ ^ ti the random variables Rit\) and Rif-i) are independent and {R(t),t > 0} 
defines a stationary zero-mean process. As usual, we denote with brackets the average taken 
over an ensemble in thermal equilibrium, i.e. ( R(t) ) = 0. 

Assuming for the frictional force the Stokes expression for a drag of spherical particle of 
radius a , we obtain the classical formula 

F v = — — V(t) , - = 6nap f u, (2.2) 
/i /i 

where \x denotes the mobility coefficient and pf and v are the density and the kinematic 
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viscosity of the fluid, respectively. If we introduce the friction characteristic time 

a — m/i , (2.3) 
the Langevin equation (2.1) explicitly reads 

^r = --V(t) + -R(t), (2.4) 
at a m 

which corresponds to the stochastic differential equation 

dV(t) = --V(t)dt + dB(t), (2.4a) 
where B(t) is a suitable chosen Wiener process. Eq. (2.4a) can be easily solved and gives 



V(t) = e- t/a [z + J e s/a dB s ^j , < t < t , 



(2.46) 



with Z = e to ^V(to). From the equation above it is easy to evaluate autocovariances; in fact, 
provided that (Z) = 0, one has, for any r, s > t : 

C v (t, s) = (V(r)V(s)) = e - (T+s)/,T | (Z 2 ) - ^ (e- 2to/a - e 2mi "(^)/-) j ) ( 2 .4c) 
where 7 is such that (B(t)B(s)} = 7 2 min(t, s), for any t,s > 0. 

At this point one usually assumes that the Brownian particle has been kept for a sufficiently 
long time in the fluid at (absolute) temperature T , so that the thermal equilibrium is reached. 
Thus, for any r, s in which the thermal equilibrium is maintained, Cy(r, s) reduces to: 

CV((r, S )) = ^e-| (2.4d) 

The equipartition law for the energy distribution requires that 

m{V 2 (t)} = kT, (2.5) 
where k is the Boltzmann constant, therefore 7 2 = 2KT jam. 

Eq. (2.4d) implies that for any t such that the equilibrium is reached the autocovariance 
functions Cy and Cr of the stochastic processes V(t) and R(t) can indeed be written as: 

Cv(t ,t) = (V(t )V(t + t)) = C v (t), t>0, (2.6) 

C R (t ,t) = {R(t )R(t + t)) = C R (t), t>0, (2.7) 
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i.e not depending on t . Moreover, by definition, the random force is uncorrelated to the 
particle velocity, namely 

C VR (t ,t) = (V(t )R(t + t)) =0, t>0. (2.8) 

Hereafter we shall assume t = 0. Because of Eqs. (2.4) and the previous assumptions, one 
has that for any t > 

C v {t) = ( V 2 (0) > e - */ " = — e - */* 7 . (2.9) 

m 

Moreover, up to a multiplicative constant 

C R (t) = — ( V 2 (0) ) S(t) = -kT6(t) , (2.10) 
a o 

where S(t) denotes the Dirac distribution. Starting from Eqs. (2.6)-(2.7), one can define the 
power spectra or power spectral densities Iv(w) and Ir(lo) , u G 1R , by taking the Fourier 
transforms of the respective autocovariance functions; i.e. 

/+oo -I />+oo 

C v (t)e- lujt dt, C v (t) = — I v (co)e +iUjt dco, (2.11) 
-oo 27T J-oo 

/ + OO i /» + oo 

C R (t)e- lujt dt, C R (t) = — I R (u)e +lujt du. (2.12) 
-oo 271 J-oo 

The result (2.9) shows that the velocity autocorrelation function decays exponentially in time 
with the decay constant a, while (2.12) means that the power spectrum of R(t) is to be 
white, i.e. independent on frequency, resulting 

I R (u) = I R = -kT. (2.13) 

(7 

It can be readily shown that the mean squared displacement of a particle starting at the 
origin at to = (displacement variance) is given by 

(X 2 (t)} = 2 f(t-T)C v (r)dT = 2 [ dn I ' dr , t>0. (2.14) 

Jo Jo Jo 

For this it is sufficient to recall that X(t) = J* * V(t') dt' , and use (2.4d). we obtain 

(X 2 (t)} = 2D t - a (l - e - ^) , *>0, (2.15) 

where 



POO 

D = a ( V 2 (0) ) = / C v (t)dt. (2.16) 
Jo 
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We note from (2.15) that 

( X 2 {t) } = 2Dt [1 - (t/a)- 1 + EST] , as t^oo, (2.17) 
(EST = exponentially small terms) so that 

D=\im^^l. (2.18) 

t— >oo 2,t 

Furthermore, using (2.3), (2.5) and (2.16), we recognize that 

D = —kT = pkT. (2.19) 

m 

The constant D is known as the diffusion coefficient and the equation (2.19) as the Einstein 
relation. 

3 The hydrodynamic approach to the Brownian Mo- 
tion 

On the basis of hydrodynamics, the Langevin equation (2.4) is not all correct, since it ignores 
the effects of the added mass and retarded viscous force, which are due to the acceleration 
of the particle, as pointed out by several authors. 

The added mass effect requires to substitute the mass of the particle with the so-called 
effective mass, m e = m [1 + Pf/(2p p )] , where p p denotes the density of the particle, see e.g. 
Batchelor (1967) [3]. Keeping unmodified the Stokes drag law, the relaxation time changes 
from a = m p to a e = m e p : thus 

<7 e = <7 (l + -L), With X = — • (3-1) 

The corresponding Langevin equation has the form as (2.4), by replacing m with m e and a 
with o~ e . With respect to the classical analysis, it turns out that the added mass effect, if it 
were present alone, would be only to lengthen the time scale (a e > a ) in the exponentials 
entering the basic formulas (2.11) and (2.15) and to decrease the velocity variance (^ 2 (0)) , 
consistently with the energy equipartition law (2.5) at the same temperature. Consequently, 
using (3.1), the diffusion coefficient is unmodified and turns out to be 

D = a e (V 2 (0)) =pkT } (3.2) 

so the Einstein relation (2.19) still holds. 

The retarded viscous force effect is due to an additional term to the Stokes drag, which is 
related to the history of the particle acceleration. This additional drag force, proposed by 
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Basset and Boussinesq in earlier times, and nowadays referred to as the Basset history force, 
see e.g. Maxey & Riley (1983) [52J, reads (in our notation) 



fi V 7T 



t>. 



dV(r)/di 



■dr, t 



— , t > i* > — oo. 
v 



(3.3) 



We suggest to interpret the Basset force in the framework of the fractional calculus. In this 
respect, taking t* = , we write 



dV(r)/dr 



dr 



d 1 ' 2 
dt 1 / 2 



V{t) 



(3.4) 



where d 1 ' 2 /dt 1 / 2 denotes the fractional derivative of order 1/2 (in the Caputo sense), see for 
details Caputo (1967, 1969) |BIE], Caputo & Mainardi (1971) [ID], Mainardi (1996, 1997) 
[Ml E], Gorenflo & Mainardi (1997) [32], and Podlubny (1999) [61]. This definition of 
fractional derivative differs from the standard one (in the Riemann-Liouville sense) available 
in classical textbooks on fractional calculus, see e.g. Oldham & Spanier (1974) [M] Ross 
(1975) [69], Samko et al. (1993) [72] and Miller & Ross (1993) p]. In fact, if /(*) denotes 
a causal (sufficiently well-behaved) function and < a < 1 we have 



) 

RL 



\dt a 



1 



dr 



r(l -a)dt I /(T) (t - r) c 



[/(*)] 



r 



(3.5) 
(3.6) 



c r(l — a) 7 dr (t 

where T denotes the Gamma function, and the suffices RL and C refer to Riemann-Liouville 
and to Caputo, respectively. 

We'recall the relationship between the two fractional derivatives (for a sufficiently well- 
behaved function): 



Oo 1/(4)1 



d° 



dt c 



RL 



f(Q + )t- a 
r(i-a) 



(3.7) 



Then, using (2.4), (3.1) and (3.3-4), the Langevin equation turns out to be 



dV 
~dt 



1 

0~e 



d 1 ' 2 
dt 1 / 2 



V(t) + —R(t) 



(3.J 
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where the suffix C is understood in the fractional derivative. We agree to refer to (3.8) as 
to the fractional Langevin equation. 

It is worth noting that if the process is meant to be in thermodynamic equilibrium (at 
to = 0), we should account for the hydrodynamic interaction long memory, and thus it is 
correct to integrate the Langevin equation (3.8) from t* = — oo . Basing on the observations 
by Dufty (1974) [21] and Felderhof (1978) [26], we introduce in (3.8) the random force 

fl,(f) = R(t) f ^pHi dr . (3.9) 

In view of the fluctuation- dissipation theorem, Kubo (1966) [13] considered a generalized 
Langevin equation (GLE), introducing a memory function j(t) to represent a generic 
retarded effect for the friction force. In our case Kubo's GLE reads (in our notation) 

dV C l 1 

4r = - l(t~r) V(t) dr + — R.(t) , (3.10) 
at Jo m e 

where R,(t) = R(t) - m e /"^ j(t - r) V(r) dr . 

For the sake of convenience, from now on we shall drop the suffix * in the Langevin equations. 
Basing on the fundamental hypothesis (2.8). i.e. ( V(0) R(t) ) = 0, t > 0, and using the 
Laplace transform, 

POO 

f{t) + f{s)= / e~ st f(t)dt, seV 



(where the sign -r denotes the juxtaposition of a function depending on t with its 
Laplace transform depending on s ), the fluctuation- dissipation theorem is readily expressed, 
according to Mainardi & Pironi (1996) [5U], as 

■ (V 2 (0)) 

C v (s) = (V(0)V(t))= { —i J U-, (3.11) 

C R (s) = ( R(0)R(t) ) = m 2 e ( V 2 (0) ) 7( a ) . (3.12) 
For the proof of (3.11)-(3.12) see Appendix A. 

We understand that the causal approach to Kubo's fluctuation-dissipation theorem can be 
controversial. For a criticism to the causal approach we refer the reader to a paper by 
Felderhof (1978) [26]. 

The classical results are easily recovered for t > noting that, in the absence of added mass 
and history effects, we get 7(5) = 1/er -r j(t) = 5(t)/a . 
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For our fractional Langevin equation (3.8), we note that 

S(t) 



700 = - [l + v^^ 1/2 ] -7W = - 
(j P dp 



+-3/2 
Zy/TV 



(3.13) 



where Q(t) is the Heaviside step function. Therefore the expression for ^(t) turns out to 
be defined only in the sense of distributions. Specifically, 8(t) is the well-known Dirac delta 
function and t~ 3 ^ 2 0(i) is the linear functional over test functions, <p(t) , such that 



^-3/2 



Q(t) , <f>(t)) 



W) - 0(0)] 

£3/2 



For more details on distributions, see e.g. Gel'fand & Shilov (1964) [30] or Zemanian (1965) 



The significant change with respect to the classical case results from the t _3//2 term. Not only 
does it imply a non-instantaneous relationship between the force and the velocity, but also 
it is a slowly decreasing function so that the force is effectively related to the velocity over a 
large time interval. The representation of the force in terms of distributions, as required by 
the GLE, is not strictly necessary since we can use the equivalent fractional form. 

Let us now consider the correlation for the random force. The inversion of the Laplace 
transform Cr(s) yields, using Eqs. (3.12)-(3.13), 



mz 



C R (t) = ^(V'(0)) 



S(t) 



f-3/2 

r — 9(t) 



7T 



(3.14) 



Thus, from the comparison with the classical result (2.12), we recognize that, in the presence 
of the history force, the random force cannot be longer represented uniquely by a white 
noise; an additional "fractional" or "coloured" noise is present due to the term £~ 3 / 2 which, 
as formerly noted by Case (1971) [T2], is to be interpreted in the generalized sense of 
distributions. 



Let us now consider the velocity autocorrelation. Inserting (3.13) in (3.11), it turns out 



^ 2 (0)> 



CviS) -—[l + V^os^]/ 
where, because of (2.2-3) and (3.1), (3.3), 

9 



V\0)) 



Op 



(Jp s 



1/2 



1/Ve 



(3.15) 







7£ 



= % 

2x + l 2 Pp + Pf 



(3.16) 



We note from (3.16) that < < 9, the limiting cases occurring for \ = oo and \ = 0, 
respectively. We also recognize that the effect of the Basset force is expected to be negligible 
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for (3 — > , i.e. for particles which are sufficiently heavy with respect to the fluid (p p ^> pf). 
In this case we can assume the validity of the classical result (2.11). 

Applying in (3.15) the asymptotic theorem for Laplace transform as s — > , see e.g. Doetsch 
(1974), HH] we get as t -»■ oo , 

C v (t) ~ ( V 2 (0) > VWN (t/a e y 3/2 , t->oo. (3.17) 

The presence of such a long-time tail was observed by [2] in a computer simulation of velocity 
correlation functions. Furthermore, from (3.15) it is easy to obtain the following results 

POO 

CV(0) = lim sC v (s) = (V 2 (0)), C v (t)dt = C v (0)=a e (V 2 (0)). (3.18) 



The explicit inversion of the Laplace transform in (3.15) can be obtained basing on Appendix 
B, see also Mainardi et al. (1995) [5TJ, Mainardi (1997) [47], and reads 



' a+E 1/2 (a + y/i) - a_£ 1/2 (a_ VI) -yffl ± {fi - A)W 



C v {t) 



a± 



a + — a_ 2J1J, 



E 1/2 (ay/i) [1 + 2a 2 t] + la^ftfix a = — — (/3 = 4), 



(V 2 (0)) 
where 

oo 

r(n/2 



(3.19) 



00 a n -j.n/2 2 

E 1/2 (a v / t) = ^ =e a hM-aVt) (3.20) 



n=0 

denotes the Mittag-Leffler function of order 1/2 and erfc the complementary error function. 
For properties of the Mittag-Leffler function we refer the reader to Erdelyi (1955) [23J, 
Gorenflo & Mainardi (1997) [32], Podlubny (1999) [M] and Mainardi & Gorenflo (2000) ggj. 

Let us now consider the displacement variance, which is provided by (2.14). From the 

Laplace transform (X 2 (s) ) = 2Cy(s)/s 2 , we derive the asymptotic behaviour of (X 2 (t) ) 
as t — > oo , which reads 

(X 2 (t)) = 2Dt{l-2 v ^(t/a e )- 1/2 + 0[(t/a e )- 1 ]^ , t -> oo , (3.21) 

where D is the diffusion coefficient (3.2). The explicit expression of the displacement variance 
turns out to be, taking j3 ^ 4 , 



X"(t) ) = 2D <t - 2^ + < 11 ~ -^>1 - a - I 1 -^Ml ^ . (3.22) 

Thus, the displacement variance is proved to maintain, for sufficiently long times, the linear 
behaviour which is typical of normal diffusion (with the same diffusion coefficient of the 
classical case). However, the Basset history force, which is responsible of the algebraic decay 
of the velocity correlation function, induces a retarding effect (oc t 1 ' 2 ) in the establishing of 
the linear behaviour. 
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4 Numerical results and discussion 

In order to get a physical insight of the effect of the Basset history force we exhibit some 
plots concerning the velocity autocorrelation (3.19) and the displacement variance (3.22), for 
some values of the characteristic parameter x — Pp/Pf ■ As we have noted before, the Basset 
retarding effect is more relevant when the parameter j3 introduced in (3.16) is big enough, 
namely when x '■= Ppl Pf is sufficiently small. In this section, devoted to numerical results 
and discussion, the plots will clearly show the increasing effect as x becomes smaller and 
smaller. But in order to recognize the type of anomalous diffusion induced by the behaviour 
of the displacement variance what can we do? We postpone our approach to this problem 
to the final discussion. 

We now agree to take non-dimensional quantities, by scaling the time with the decay constant 
a of the classical Brownian motion and the displacement with the diffusive scale (Da) 1 / 2 . 
With these scales the asymptotic equation for the displacement variance reads ( X 2 (t) ) ~ 2 t . 

For decreasing values of x — 1,0.5,0.1 we consider the velocity autocorrelation normalized 
with its initial value ( V 2 (0) ) and the displacement variance normalized with its asymptotic 
value 2 1 . In Figure 1 we compare versus time the functions Cy and ( X 2 ) / (2t) provided by 
our full hydrodynamic approach (added mass and Basset force), in continuous line, with the 
corresponding ones, provided by the classical analysis, in dashed line, and by the only effect 
of the added mass, in dotted line. For large times we also exhibit the asymptotic estimations 
(3.17) and (3.21), in dotted line, in order to recognize their range of validity. 

The correlation plots exhibit the well-known algebraic tail. The time necessary to reach 
the asymptotic behaviour increases as the density ratio x decreases. By comparing the two 
figures, it appears that the variance approaches the asymptotic regime as the autocorrelation 
becomes sufficiently small, independently on its time dependence. We note that in the time 
interval necessary to reach the asymptotic behaviour the displacement variance exhibits a 
marked deviation from the standard diffusion. Because this time interval turns out to be 
orders of magnitude longer than the classical one, it is relevant to discuss about various 
diffusion regimes before the normal one is established and consequently about the nature of 
the anomalous diffusion. 

In order to explore the nature of the anomalous diffusion we follow the following reasonings. 
At first, we find it convenient to introduce the instantaneous diffusion coefficient 

We note that such coefficient is suitable to be measured in experiments. 

If we assume that (X 2 (t) ) varies in a range of t, say t k < t < t k+ i, approximatively as a 
power law 

(X 2 (t)) ~at a \ t k <t<t k+1 , (4.2) 
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Figure 1: From top to bottom plots on the velocity correlation and on the displacement 
variance for x — 1, 0.5, 0.1 
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where a is assumed to be a suitable constant of dimensions L 2 T a , called effective diffusion 
coefficient, we are usually led to say that the diffusion regime is anomalous if a ^ 1, precisely 

• slow diffusion or sub-diffusion if < < 1 , 

• fast diffusion or super- diffusion if 1 < < 2 . 

When «fc = 2 the diffusion regime is said to be ballistic; and, of course, when a/. = 1 the 
diffusion is normal. 

We easily recognize from Eq. (4.2) that the power a& can be written as 

d(log(X 2 (t))) 



h<t<t 



d(logt) 

In fact from (4.2) we get for any tk <t < t^+i 

log ( X 2 (t) ) = log [a t Qfc ] = log a + a fc log t 



k+l- 



(4.3) 



(4.3) 



Then, varying k we are able to define a function a(t), representing a sort of instantaneous 
diffusion parameter, which expresses, locally, the characterization of the regimes of 
anomalous diffusion. In the following plots we shall exhibit (X 2 (t)) as given by (3.22) 
with D(t) as given by (4.1) and in addition log (X 2 (t) ) with a(t) as given by (4.3). 
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We will note that in all cases D(t) is an increasing function of time, starting from at t = 
and getting 1 in the limit as t — > oo. On the other hand, a(t) is a decreasing function of 
time, starting from 2 at t = (ballistic regime) and getting 1 in the limit as t — > oo (normal 
diffusion) . 

In our plots we scale the time t with the equilibrium relaxation a e given by (3.1). and, as 
usual, we assume for \ lower and lower values x — 1,0.5,0.1. We easily recognize that the 
retarding effect of the Basset force is larger for lighter particles, so for them the exponent 
a(t) remains over 1 for longer times so that, in experiments of limited time duration, this 
effect could be interpreted as fast diffusion. 

5 Conclusions 

The velocity autocorrelation and the displacement variance of a Brownian particle moving in 
an incompressible viscous fluid are known to be the fundamental quantities that characterize 
the evolution in time of the related stochastic process. Here they are calculated taking into 
account the effects of added mass and both Stokes and Basset hydrodynamic forces that 
describe the friction effects, respectively in the steady state and in the transient state of the 
motion. 

The explicit expressions of these quantities versus time are computed and compared with 
the respective ones for the classical Brownian motion. 

The effect of added mass is only to modify the time scale, that is the characteristic relaxation 
time induced by the Stokes force. The effect of the Basset force, which is of hereditary type 
namely history-dependent, is to perturb the white noise of the random force and change the 
decay character of the velocity autocorrelation function from pure exponential to power law 
because of the presence of functions of Mittag-Leffter type of order 1/2. 

Furthermore, the displacement variance is shown to exhibit, for sufficiently long times, the 
linear behaviour which is typical of normal diffusion, with the same diffusion coefficient of 
the classical case. However, due to the Mittag-Leffler functions, the Basset history force 
induces a very long retarding effect in the establishing of the linear behaviour, which could 
appear, at least for light particles, as a manifestation of anomalous diffusion of the fast type 
(super- diffusion) . 
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Appendix A 

Let us consider the generalized Langevin equation (3.10), that we write as 

R(t) = m e \v(t) + 7 (t) * V(t) 



(ai; 



where • denotes time differentiation and * time convolution. The assumption of stationarity 
for the stochastic processes along with the following hypothesis 

(R(t)) = 0, ( V(0) R(t) ) = , t>0, (A2) 

allows us to derive, by using the Laplace transforms, the two fluctuation- dissipation theorems 

(V 2 (0)} 



C v (s) := (V(O)V(t)) 



s + 7(s) 



(A3) 



and 



C R (s) :=(R(0)R(t)) = ml(V 2 (0))j(s). 



(A4) 

Our derivation is alternative to the original one by Kubo (1966) [43J who used Fourier 
transforms; furthermore, it appears useful for the treatment of our fractional Langevin 
equation. 

Multiplying both sides of (A.l) by V(0) and averaging, we obtain 



(V(0)V(t))+y(t)*(V(0)V(t)) = 0. 
The application of the Laplace transform to both sides of (A. 5) yields 

s ( v(oyv(t) > - ( v 2 (o) > + 7( a ) ( v(oyv(t) ) = o , 

from which we just obtain (A. 3). 

Multiplying both sides of (A.l) by R{0) and averaging, we obtain 

C R (t) := ( R(0) R(t) ) = m\ [( V(0) V(t) ) + 7 (t) * ( V(0) V(t) ) 

Multiplying both sides of (A.l) by R(0) and averaging, we obtain 

C R (t) := ( R(0) R(t) )=ml[( V(0) V{t) ) + 7 (t) * ( V(0) V(t) ) 

Noting that, by the stationary condition, 

( V(0) V(0) ) = , ( V(0) V(t) } = -( V(0) V(t) ) , 



(A5) 



(A6) 



(A7) 



(A7) 



(A8) 
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the application of the Laplace transform to both sides of (A. 7) yields 

C R {s) = ml {s(V^)V{t))-l{s) [s(V{0)V(t))-(V\0) 

Since 



( 1/(0) V(t) ) = -( V(0) V(t) ) = -s( V(0) V(t) ) + ( V 2 (0) ) , 



we get 



(A9) 
(A10) 
(All) 



C R (s) = m 2 e {s [s C v {s) + ( V 2 (0) ) - 7( a ) C v (s) ] + j(s) ( V 2 (0) ) } , 
from which, accounting for (A. 3), we just obtain (A. 4). 

Appendix B 

In this Appendix we report the detailed manipulations necessary to obtain the result (3.19) 
as Laplace inversion of (3.15). For this purpose we need to consider the Laplace transform 



N(s) 



1 



s + bsW + l ' 



and recognize that 

C v (t) 



<V*(0) 



= N(t/a e ) -T- o e N(a e s) = 



(5.2) 



where we have used the sign for the juxtaposition of a function depending on t with 
its Laplace transform depending on s . The required result is obtained by expanding N(s) 
into partial fractions and then inverting. Considering the two roots A± of the polynomial 
P{z) = z 2 + bz + 1 with z = s 1 / 2 , we must treat separately the following two cases: 
i) < b < 2 , or 2 < 5 < 3 , and n) 6 = 2, which correspond to two distinct roots 
(A+ 7^ A_), or two coincident roots (A + = A_ = —1), respectively. We obtain 



b^2 



/3^4, ^5/8, 
N(s) = 1 



A 



A. 



S + b S V2 + 1 s l/2 ( s l/2 _ A+ ) s l/2 ( s l/2 _ A _) 



with 



ii) b 



A 



± 



-6 ± (b 2 - 4) 



1/2 



A- 



A, 



A i - A_ 



(3 = 4, X = 5/8, 
iV( S ) 



S + 2 S V2 + 1 ( S V2 + 1)2 



(B.3) 



(BA) 



(B.5) 
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The Laplace inversion of (B.3) and (B.5) turns out, see below, 



/ ) A + E 1/2 (X + Vt) - A_E 1/2 (\_Vi), 



where 

00 \n f n/2 2 

E 1/2 (X Vt) = g r(w/2 - 1} = e A * erfc(-A v^) (5.7) 

denotes the Mittag-Leffler function of order 1/2 and erfc denotes the complementary error 
function. In view of (B.l-2), equation (B.6) is equivalent to (3.19). 

Let us first recall the essentials of the generic Mittag-Leffler function in the framework of 
the Laplace transforms. The Mittag-Leffler function E a (z) with a > , so named from the 
great Swedish mathematician who introduced it at the beginning of this century, is defined 
by the following series representation, valid in the whole complex plane, 

E a (z)=J2 r( Z — 7T , «>0, ze(D. (B.8) 
' 1 (an + 1) 



n=0 



It turns out that E a (z) is an entire function, of order p = 1/a and type 1 , which provides 
a generalization of the exponential function. 

The Mittag-Leffler function is connected to the Laplace integral through the equation 

,00 1 

e~ u E a (u a z) du= - , a>0. (B.9) 

Jo 1 — z 

This integral is fundamental in the evaluation of the Laplace transform of E a (A t a ) with 
A G C and t > 0. Putting in (B.9) u = st and u a z = At", we get the following Laplace 
transform pair 

E a (\t a ) + * _ - , J Res>|A| 1/a . (5.10) 

We note that, up to our knowledge, in the handbooks containing tables for the Laplace 
transforms, the Mittag-Leffler function is ignored so that the transform pair (B.10) does not 
appear if not in the special case a = 1/2 . In fact, in this case we recover from (B.10) the 
basic Laplace transform pair 

^(../.-^ ^M), (AH) 

where the Mittag-Leffler function can be expressed in terms of known functions, as shown 
in (B.7). As an exercise we can derive from (B.ll) the following transform pairs and 
consequently the result (B.6): 

' ! + ,^<L u + -L + A£ 1/2 (Av^), (B.12) 



S 1/2_ A s l/2 S 1/2( S 1/2_ A ) • 
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1 1 A 

( S V2 - A) 2 ~ S l/2 ( s l/2 _ A ) + s l/2 ( s l/2 _ A )2 

- 2Ay / ^+(l + 2A 2 t)E 1/2 (A v / t) 
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